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Cubic  Spline  Approximation  Techniques 
for 

Parameter  Estimation  in  Distributed  Systems 

H.  T.  Banks,  .J.  M.  Crowley,  and  K.  Kunisch 

ABSTRACT 

Approximation  schemes  employing  cubic  splines  in  the  context 
of  a  linear  semigroup  framework  are  developed  for  both  parabolic  and 
and  hyperbolic  second  order  partial  differential  equation  parameter 
estimation  problems.  Convergence  results  are  established  for  problems 
with  linear  and  nonlinear  systems  and  a  summary  of  numerical  experiments 
with  the  techniques  proposed  is  given. 


1.  Introduction. 


In  many  problems  of  practical  importance  one  is  confronted 
with  the  task  of  estimating  unknown  parameters  in  mathematical  models 
from  certain  observations  of  the  underlying  physical,  biological,  etc. 
phenomenon  being  modeled.  This  paper  is  devoted  to  the  study  of  spline 
approximation  techniques  for  the  identification  or  estimation  of  constant 
parameters  in  partial  differential  equations.  The  results  presented 
here  are  actually  the  outcome  of  investigations  partly  reported  in  [9] 
in  which  a  theoretical  convergence  framework  was  developed  and  applied 
to  treat  "modal"  approximation  schemes  for  identification  and  control. 
Here  we  use  the  theoretical  framework  of  [9]  to  treat  spline  based 
techniques  that  were  developed  and  tested  simultaneously  with  the  "modal" 
methods  of  [9] . 

Two  specific  classes  of  problems  are  investigated  below.  In 
section  2  we  treat  an  identification  problem  for  a  class  of  parabolic 
equations,  whereas  results  for  hyperbolic  equations  are  given  in  section 
4.  Sections  3  and  5  are  devoted  to  discussions  of  the  implementation 
of  the  approximation  schemes  along  with  numerical  examples  involving 
parabolic  and  hyperbolic  equations  respectively.  In  both  cases 
interesting  technical  questions  arise.  If  the  parabolic  equations 
arc  studied  as  a  system  in  the  usual  state  spaces,  then  fit-to- 

data  criteria  involving  only  integral  terms  can  be  treated  with  relative 
ease,  whereas  fit-to-data  criteria  that  employ  point  spatial  evaluations 
(which  are  often  used  in  practice)  present  some  essential  difficulties. 

In  the  case  of  hyperbolic  equations,  the  semigroup-theoretic  approach 


that  is  taken  in  this  paper  leads  to  a  state  space  with  inner  product 
depending  on  one  of  the  unknown  parameters  of  the  equation.  For  this 
reason  it  is  desirable  (in  other  problems  this  is  essential  -  see  the 
discussions  regarding  certain  unknown  function  space  parameters  in 
[4]  and  [9])  to  develop  parameter  dependent  convergence  results  in  the 
context  of  a  parameter  dependent  state  space  framework. 

In  our  presentation  below  a  number  of  numerical  examples  are 
given  to  illustrate  the  theoretical  convergence  results.  For  a  survey 
of  the  potential  uses  of  our  semigroup- theory  based  parameter  estimation 
techniques  to  problems  in  a  number  of  specific  areas  of  applications 
(e.g.,  reservoir  and  seismic  engineering,  large  space  structures, 
transport  models  in  physiology  and  population  biology,  elasticity  and 
others)  the  interested  reader  is  referred  to  [4]  and  a  monograph  that 
the  authors  currently  have  in  preparation. 

We  defer  until  section  6  below  a  number  of  comments  on  other 
work  related  to  the  ideas  presented  here. 

The  notation  used  throughout  this  paper  is  rather  standard  and 
follows  closely  that  explained  in  [9].  For  norms  of  elements  in  Banach 
spaces  we  use  |*|, whereas  ||  »  |[  is  used  for  operator  norms.  We  shall 
occasionally  use  a  subscript  such  as  ]*|g  to  distinguish  a  norm  in  the 
space  B  from  other  norms.  The  standard  practice  of  denoting  the  usual 
norm  in  Lp(0,l),  p=2,°°,  by  |»|  will  be  followed. 


-S- 


2.  A  class  of  parabolic  partial  differential  equations. 

As  a  first  example  for  the  techniques  that  were  mentioned  in  the 
introduction  we  consider  the  heat  equation 

q1 

(2.1)  ut  =  (y-)D(p  Du)  +  q2u  +  f(q3,t,x,u) ,  for  t  >  0,  x  €  [0,1]  , 
with  boundary  and  initial  conditions, 


(BC) 

u(t ,0)  =  u ( t , 1 )  =  0, 

for  t  > 

0, 

(IC) 

u(0,x)  =  q4<p(x) ,  for 

x  6  [0, 

1], 

where 

u  =  u(t,x)  €  R  and 

II 

a 

We  assume  that  (2 . 1) - (BC) - (IC) 

models  a  phenomenon  for  which  we  wish  to  "identify"  or  estimate  the 

parameter  vector  q  =  (q^,...,q^)  £  Q  cR^  from  known  measurements 

?  =  {y(t.))r  ,  with  y(t.)€  Rm  taken  at  times  t.  in  the  fixed  interval 
1  i=l  1  v  l  i 

[0,T]  of  observation.  These  observations  y(t^)  correspond  to  values 
C C t ^ ,q)£  (t^ ,q)  in  the  mathematical  model  where  C(t^,q)  is  a  real 

m  * l  matrix  which  is  continuous  in  q  for  each  t^,  and 

C(tt,q)  =  col(u(t. ,x  ;q) , . . . ,u(t. ,x^;q)) . 

Subsequently  we  shall  give  conditions  on  f  that  guarantee  the  existence 
of  a  solution  u  of  (2. 1) - (BC) - (IC) ,  which  will  also  be  denoted  by 
u(*,*;q)  whenever  we  wish  to  emphasize  dependence  on  q.  Since  the 
problem  of  determining  a  vector  parameter  q  €  Q  such  that 
y(t^)  =  C (t^ ,q)4 (t^ ,q)  for  all  i  would  most  often  lead  to  an  unreasonable 
and  mathematically  ill  posed  problem,  we  formulate  the  problem  of  finding 
an  estimate  q  as  an  optimization  problem  in  the  following  way: 

r  2 

(ID)  Minimize  J(q,u(* , • ;q) ,y)  =  £w  | y (t . )  -  C(t. ,q)5 (t. ,q) | 

i=l  1  1 
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4 

over  q  £  Qc  R  subject  to  u(*  ,*  ;q)  satisfying 
(2.1)-(BC)-(IC)  on  [0,T]  *  [0,1] 

In  the  above  formulation  of  the  cost  functional  the  weights 

w.  >  0  can  be  used  if  needed  to  compensate  for  a  priori  known  extreme 

behavior  of  the  solution,  as  for  example,  exponential  growth  or  decay. 

(.For  simplicity,  we  take  w.  =  1  in  all  our  discussions  and  examples  in  this  paper 
The  cost  functional  that  is  chosen  here  is  but  one  of  the  several  that 
are  of  practical  relevance.  The  one  that  we  chose  for  our  presentation, 

however,  exhibits  many  interesting  technical  difficulties  (due  to  the 

use  of  point  evaluations);  in  fact,  our  final  convergence  result,  when 

used  with  this  functional,  requires  restrictions  on  f  (dependence  on 

t,  x  is  allowed,  but  no  dependence  on  u  is  permitted).  For  the  case 

of  cost  functionals  involving  distributed  measurements  as, for  example, 
r  fl  ^ 

J(q,u,y)  =  I  | u (t . , x ; q)  -  y(t.,x)|  dx,  our  approximation  results  are 
i  =  lJn  1  1 

valid  for  equations  containing  quite  general  nonlinearities  f  (see 

[  9]). 

Before  we  formulate  the  identification  problem  in  the  Hilbert 
space  X  =  H^  with  inner  product  f  k(x)q>(x)i(»(x)dx  and  asso- 

dated  norm  [•[,  we  summarize  the  hypotheses  that  will  be  needed  through¬ 
out  this  section.  We  denote  by  H1  the  usual  Sobolev  spaces  over  (0,1) 
as  discussed  in  [  1  ].  The  function  F:  Q  *  [ 0 ,«>)  *  H°  ->  H°  which  will 
be  used  below  is  defined  to  be  the  composition  map  F'q,t,v)  =  f(q,t, *,v(»)) • 
Further,  for  a  given  parameter  set  Q,  we  let  Q4  =  (q^:  q  =  (qj,...,q^)  £  Q}. 

Our  hypotheses  are  listed  for  easy  reference  as  follows: 

(HQ)  Q  is  a  compact  subset  of  R4  and  there  exist  positive 

numbers  q^  and  q^  such  that  q  €  Q  c  R4  implies 
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q,  i  q,  £  q 


1  • 


(HP)  The  functions  p  and  k  satisfy  p  €  C  [0,1],  k  €  C[0,1], 
with  k(x)  >  0  and  p(x)  >  0  for  all  x  €  [0,1]. 

(HF)  The  nonlinear  function  F  satisfies 

(i)  for  each  continuous  function  u:[0,T]  -+•  X  and  each 
q  €  Q,  the  map  t  -v  F(q,t,u(t))  is  measureable, 

(ii)  for  each  constant  M  >  0  there  exists  a  function 
^1  =  such  that  for  any  q  €  Q  we 

have 

(FCq.t.Uj)  -  F(q,t,u2)|  <_  kj  (t)  ]  Uj-u^ 
for  all  u1«u2  €  X  with  |u.|  £  M, 

(iii)  there  exists  a  function  k2  in  L2(0,T)  such  that 
|F(q,t,v)|  k2(t){|v|  +  1), 

for  all  v  €  X,  q  6  Q  , 

(iv)  for  each  (t,v)  €  [0,T]  *  X  the  map  q  -+•  F(q,t,v) 
is  continuous. 

In  [  9 ]  we  have  detailed  conditions  on  the  perturbation  f 
that  guarantee  (HF) .  For  our  parabolic  systems,  statement  of  the  con¬ 
ditions  of  (HF)  in  terms  of  f  are  rather  obvious  except  in  the  case 
of  (HF)(ii)  which  requires  a  global  Lipschitz  criterion  for  f  with 
respect  to  u  due  to  the  fact  that  we  are  using  the  l2  norm  in  the 
local  statement  of  (HF)(ii)  for  F.  We  note  that  as  a  consequence  of 
(HP)  the  spaces  X  and  H^  have  equivalent  norms. 
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We  shall  consider  the  semi  linear  equation 

=  A(q)u(t)  +  F (q , t , u(t) )  for  t  >  0, 

(2.2)  UT 

u(0)  =  uQ (q)  =  q , 

where  ACq)^  =  O(pD^)  +  q^ ,  with  dom  A(q)  =  {^  €  H  (0)  =  ip(l)  =  0}. 
Under  (HQ)  it  is  well  known  that  A(q)  is  the  infinitesimal  generator 
of  a  linear  ^-semigroup  T(t;q)  for  each  q  £  Q.  We  point  out  that 
we  have  used  the  letter  u  in  two  different  ways:  In  (2.1)  u  =  u(t,x)eR, 

whereas  in  (2.2)  u(t)  or  u(t;q),  if  dependence  on  q  is  emphasized,  is 
an  element  of  X.  This  should  not  create  any  ambiguities,  since  the 
notion  of  solutions  of  (2. 1)- (BC) - (IC)  or  (2.2)  will  be  fixed  throughout 
the  paper  to  be  that  of  mild  solutions.  We  recall  that  t  -►  u(t;q)  is 
called  a  mild  solution  of  (2.2)  on  [0,T]  if  it  satisfies 

(2.3)  u(t;q)  =  T(t;q)uQ(q)  +  T(t-s;q)F(q,s,u(s;q))ds, 

b 

for  t  £  [0,T].  Under  hypothesis  (HF)  one  can  easily  demonstrate 
existence  of  a  unique  solution  to  (2.3Xe.g.  see  [9]).  The  relationship 

between  mild  and  strong  solutions  of  (2.2)  and  classical  solutions  of 
(2 . 1 ) - (BC) - ( IC)  has  been  the  focus  of  many  investigations  and  we  onlv  refer 
to  [15]  as  one  possible  reference.  We  are  now  prepared  to  formulate  the 
abstract  identification  problem  associated  with  equation  (2.2): 

(IDA)  Minimize  J (q,u(. ;q) ,y)  over  q  £  Q  c  R4  subject  to 
u(.;q)  satisfying  (2.3)  on  [0,T]. 

An  unavoidable  difficulty  with  parabolic  problems  now  becomes 
apparent.  We  consider  (2.2)  as  an  equation  in  ii^  and  a  cost  functional 
J  which  involves  point  evaluation,  an  operation  that  is  not  well-defined 
on  H  .  Moreover,  it  will  shortly  become  evident  that  "qN  -+  q*  implies 
u  (t;q  )  -►  u(t;q*)",  where  u  (t;q)  is  an  approximation  to  u(t;q)  in 


-9- 


X,  is  an  analytical  statement  of  central  importance  in  our  approach  to 

N  N 

the  parameter  identification  problem.  Convergence  of  u  (t;q  )  to 

u(t;q*)  in  X,  however,  clearly  does  not  imply  the  desired  convergence 

N  N  N  ^ 

of  the  cost  functionals  J(q‘,u  (• ;q  ) ,y)  to  J(q* ,u(* ;q*) ,y) .  As  a 

quick  remedy  to  these  difficulties,  one  might  be  tempted  (this  is  not 

our  remedy)  to  set  up  (2 . 1) - (BC) - (IC)  in  a  state  space  with  a  stronger 

topology.  This  would  involve  a  more  complicated  inner  product,  of 

course,  and  consequently  lead  to  a  more  complex  structure  for  the  matrix 

representations  of  the  approximating  finite  dimensional  problems,  a 

feature  that  is  highly  undesirable  from  the  numerical  point  of  view. 

With  regard  to  the  point  evaluations  used  in  J ,  we  note  that 

as  a  consequence  of  the  "smoothing  effect"  of  parabolic  systems,  under 

(HF)  one  can  use  the  theory  of  monotone  operators  to  argue  that 

1/2  1 

u(t;q)  6  dom  (u)I  -  A(q))  c  H  ,  for  some  >  0.  Thus  point  evaluation 
can  be  justified;  this  is  discussed  in  more  detail  in  [9]-  For  a 
solution  to  the  second  of  the  above-mentioned  difficulties  see  Theorem 
2.1  below. 

We  now  explain  a  Galerkin  approach  employing  cubic  spline 
subspaces  to  solve  (IDA)  iteratively.  Given  a  value  of  N  and  a  vector 
q,  we  seek  an  approximate  solution  to  (2.2)  in  X^  =  span(Bp, . . . ,B^} 
of  the  form 

(2.4)  uN(t;q)  =  l  w"(t)B*  =  Z  w"(t;q)B"  , 

N  j=0  J  J  j=0  J  J 

where  (Bj)  is  the  set  of  cubic  spline  basis  functions  appropriately 

modified  to  be  in  dom  A(q).  More  precisely,  let  AN  =  with 

i  ~N 

xi  =  N  ’  ^or  1  =  am*  let  Bj >  J  =  -1 . N+l ,  denote  the 
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standard  Cr  (0,1)  basis  elements  for  the  cubic  B-spline  subspaces 

of  dimension  N+3  with  respect  to  the  grid  A^  (see  p.  208-209  of  [17]). 

Then  is  given  by 

J  N  ~N 

B  =  B  ,  for  2  <  j  <  N-2  , 

J  J  -  J  “ 

N  ~N  ~N  N  ~N  ~N 

B0  ■  Bo  -  4B-1'  BS  '  bn  -  4C  • 

N  ~N  ~N  N  ~N  ~N 

B1  =  Bo  -  4B1  ■  bn-i  =  bn  -  4Bn-i  • 

Note  then  that  XN  =  SB(AN)  =  {<*>£  S3CAN) : (0)  =  4» ( 1 )  =  0}, 

3  N  2 

where  S  (A‘ )  =  (4>  €  C  (0,1):  <P  is  a  cubic  polynomial  on  each  interval  [x.  ,x^+^]} 


(2.5) 


The  approximate  solutions  to  (2.2)  are  determined  from  requiring 
N 

that  for  all  2  €  X  ,  we  have 

<UN(t),z>=  A(q)UN(t),Z^  +  ^F(q,t,uN(t)),Z^>  , 

<^uN(0),zy=  ^q4<t>,Z^>, 

or,  equivalently, 

l  <w?(t)B?,  B*>  -  l  <w^(t)A(q)B^,  B*  >  +  <V  (q ,  t  ,UN  (t)  )  ,  B^  >  , 
j=0  J  '  1  j=o  J  • 

j=0  J  J 

This,  in  turn  can  be  written  in  matrix  form  as 

N  •  N  N  N  N  N 

Q  w  (t)  =  K  wN (t)  +  RN(F(q,t,uN(t))), 


(2.6) 

where 


N  N  N 

Q  w  (0)  =  RIN(q44>)  , 

A.,  ■  **>■ 

Aj  -  <B»,  Afq,B«> 

(R%.  =  Bj  ^  , 
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„  ,  N  .  ,  N  N  N, 

and  w  =  col(wQ,  Wj,. . . ,wN) 


Alternatively,  (2.6)  can  be  viewed  as  the  result  of  projecting 
the  original  problem  (2.2)  onto  the  finite  dimensional  subspaces  X^. 

In  terms  of  projections,  we  can  write  the  Galerkin  equations  in  the 


(2.7) 


•N  N  N  N  N 

w"(t)  =  A”(q)wJV(t)  +[P  F(q,t,u  (t))J,  for  t  >  0, 

wN(0)  =  [PN(q/)], 


N  N 

where  P  :  X  -*■  X  denotes  the  canonical  orthogonal  projection  along 

jvj  X  N  ^  m  jy 

(X'  )  (  A  (q)  stands  for  the  matrix  representation  of  A  =  P  A (q ) P 

N 

and  [•]  denotes  the  coordinate  vector  for  an  element  in  X 


Thus,  AN  =  (QN)_1KN  and  [PN4>]  =  (QiN) "V (<!>)  for  U  X. 


.Ns  -1„N . 


Since  A  (q)  is  a  bounded  operator,  it  clearly 

generates  a  semigroup  T  (t;q)  =  e  ^  .  Moreover  u  (t;q)  will 

satisfy  (2. 5) (equivalently,  (2.6))  on  [0,T]  if  and  only  if  it  satisfies 

(2.8)  uN(t;q)  =  TN(t;q)PNu  (q)  +  TN(t-s ;q) PNF (q , S ,uN(s ;q) )ds . 

u  J0 

N 

Standard  Picard  iteration  arguments  yield  that  unique  solutions  u  (»;q) 

of  (2.8)  exist  under  hypotheses  (HF) . 

In  light  of  the  above  discussions,  we  therefore  formulate  the 

approximate  identification  problems  as: 

(IDAN)  Minimize  J(q,uN(- ;q) ,y)  =  I  |y(t.)  -  ,q)SN(t .  ,q)  | 

N  i  =  1 

subject  to  u‘  (*;q)  satisfying  (2.8)  on  [0,T],  where 
5N(ti,q)  =  col(uN(ti,x1;q),...,uN(t.,xJl;q)). 

Using  the  continuous  dependence  with  respect  to  q  of  C(t^,q)  and 
u^Ct^q)  for  each  i  =  l,...,r,  one  can  easily  argue  continuity  of 
q  -*■  J(q,u^(*  ;q)  ,y)  and  consequently,  by  (HQ),  there  exists  for  each 
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-N  N 

N  a  solution  q  of  (IDA  ) .  By  compactness  of  Q  we  can  extract  a 
convergent  subsequence,  again  denoted  by  q^,  with  q^  ■*  q  €  Q.  In  the 
remainder  of  this  section  we  shall  concentrate  on  proving  that  for  any 
arbitrary  sequence  (qN)  c  Q, 

(2.9)  lim|qN  -  q*j  =  0  implies  lim|uN(t;q^)  -  uttjq*)^  =  0 

N-*»  N-*00 

or  (in  the  event  one  employs  J  for  the  fit-to-data) , 

(2.10)  lim|qN  -  q*|  =  0  implies  limj uN(t;qN)  -  u(t;q*)|x  =  0 

N-*00  N-+00 

for  each  t  £  (0,T].  We  remark  that  for  the  problems  considered  in  this 

paper  we  have  u^(t;q^)  €  and  u(t;q*)  €  so  that  the  norms  in 

(2.9)  and  (2.10)  can  be  employed  without  loss  of  meaning. 

Once  (2.9)  (respectively  (2.10))  is  verified  it  follows  immediately 

that  q  is  a  solution  of  (IDA)  (respectively  of  (IDA)  with  J  replaced 

N  N  N  ^  A 

by  J)  and  that  lim  J(q‘  ,u  (•  ;q  ),y)  =  J(q,u(* ;q) ,y)  (respectively, 

N  A 

lim  J(  qN,uN(* ;qN) ,y)  =  J  (q  ,u(* ;q),y)).  Indeed  J(qN,uN(*;qN),y)  £ 

N 

N  ^ 

J(q,u  (• ;q) ,y)  for  all  q  €  Q  and  all  N.  Under  (2.9)  we  have  that 

limjuN(t;qN)  -  uftjq)^  =  0  and  lim]uN(t;q)  -  uftjq)^  for  each  t  €  [0,T] 

N->oo  N-x» 

and  q  €  Q.  Consequently,  taking  limits  in  the  above  inequality  we  obtain 

A  A 

J (q  ,u(*  ;q)  ,y)  £  J (q ,u( • ; q),y)  for  each  q  £  Q,  so  that  q  is  a  minimizer 

/\ 

of  q  •>  J(q,u(*  ;q)  ,y)  . 

N 

Since  u(t;q)  and  u  (t;q)  are  solutions  of  (2.3)  and  (2.8) 

respectively,  (2.10)  will  follow  easily  by  the  Gronwall  lemma  and  Lcbesgue's 
bounded  convergence  theorem,  once  we  have  shown 

(Hi)  j|TN(t;q)]i  _<  MeWt  ,  with  M  and  w  independent  of  N  and  q, 

(Hii)  PN  ■*  I  strongly  in  X,  and 

(Hiii)  TN(t;qN)  -*■  T(t;q*)  strongly  in  X  and  uniformly  in  t  €  [0,T], 

N 

for  any  convergent  sequence  q  q*. 


For  details  in  a  more 


general  setting,  we  refer  to  [9  ].  Due  to  the  special  choice  of 

=  PNAPN  we  find  that  (Hi)  trivially  holds:  indeed  ||T(t;q)||  <  ewt 
with  w  =  max  q2  .  Then,  since  ^A^  =  ^AP^V,  W|*P|2  it 

immediately  follows  that  (Hi)  holds.  The  usefulness  of  this  classical 

argument  in  the  context  of  general  approximation  schemes  for  dissipative 

operators  was  apparently  first  noted  in  [  7 ] .The  strong  convergence  of 
N 

P  to  I  is  established  with  the  aid  of  three  lemmas.  The  first  is 
a  statement  of  the  standard  cubic  spline-interpolation  error  bounds,  the 
second  a  statement  of  the  useful  Schmidt  inequality  and  finally  in  the 
third  lemma  we  will  establish  (Hii)  on  a  dense  subset  of  X.  That  (Hi i ) 
obtains  then  follows  from  the  uniform  boundedness  of  the  sequence  of 
projections  {P  }  . 

4  N 

Lemma  2.1.  [2Qp.54],If  z  £  H  ,  and  I  z  denotes  the  cubic  spline  inter- 
polant  of  z  in  S^(A^),  then 

|z  -  INz|2  <  C0N'4fD4z|2  , 

|D(z  -  rNz) | 2  <  c1n_3|d4z(2 

|D2(z  -  INz)|2  <  C2N'2jn4z|2  , 
where  the  C.  are  constants  independent  of  z  and  N. 


Lemma  2.2.  [20, p. 7].  If  p^  is  a  polynomial  of  degree  n 

[a,b],  then 


=1,2,3  on 


fb 

fDP 

L  n 


(x)]  dx  <  kn(b-a) 


,(x)  dx  , 


with  k  a  constant, 
n 
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Lenuna  2.3.  If  z  €  fl  H4,  then 


PNz  -  z!  <  CqN'4|d4z|  2  , 


D(PNz  -  2)|  <  CjN--* | D^z  |  2  , 

^  ^2*’  i  ‘•i  2  * 


.*-3 .  „4_ 


D2(PNz  -  z)|  <  C  N'2|d4z 


where  the  are  constants. 

14  N  N 

Proof.  Let  z  £  Hq  n  H  and  denote  by  I  z  and  I^z  the  interpolating 

3  N  3  N 

cubic  B-spline  for  z  in  S  (A  )  and  Sq(A  )  respectively.  Using  the 
boundary  conditions  we  have  |  P^z-z)  <  jl^z-z  |  =  |l^z-z  j  < /f  CqN-4  |D4z  |2  . 

where  V  =  max{k(x)|o  £x  ^  l) .  This  implies  the  first  estimate  with 


C0  =  C0 


/max  k(x)V^^ 
\min  k(x) y  ’ 


the  extremes  being  taken  over  x  £  [0,1].  To  verify 
the  second  estimate  we  make  use  of  Lemmas  2.1  and  2.2: 


D(PNZ  -  z)  j  x  <  A  |D(PNZ  -  z)  1 2  <  /T  !d(PNz  -  ljz)l2  +  ’/T 
<  /Ik 3  N  |PNz  -  I^z  1 2  +  /T|  D(INz  -  z)|  2 


DOoZ  -  z)l2 


1  /kk^  N((pNz  -  z  1 2  +  |  z  -  Ipz|2)  +  /T|  D(INz  -  z)|  2 

1  2  /ki^  n'3c0|d4z!2  +  /Tc1n'5|d4zI2  . 

This  establishes  the  convergence  of  the  first  derivative  of  the  proiected 
elements.  The  final  estimate  can  be  argued  in  a  similar  manner. 

Having  thus  established  (Hii)  we  next  turn  to  (Hiii).  It  is 
convenient  to  use  the  following  version  of  the  Trotter-Kato  theorem. 

Proposition  2,1  [  14] .  Let  (6,|-|)  and  C^,  |  •  |N) ,  N  =  1 ,2, . . . ,  be 
Banach  spaces  and  let  8^8^  be  bounded  linear  operators.  Assume 
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N 

further  that  T(t)  and  T  (t)  are  linear  Cq  -  semigroups  on  B  and 

N  ~ 

8  with  infinitesimal  generators  A  and  A  respectively.  If 

i)  lim  I |  =  |  z|  for  all  z  £  8, 

ii)  there  exist  constants  M,^  independent  of  N  such  that 

||rN(t)  (|N  <  Mewt  ,  for  t  0, 

iii)  there  exists  a  set  Q  <z  8  such  that  p  c  dom(X),  p  =  B  and 

(XQ  -  =  8  for  some  X^  >  0  and  for  which  for  all  z  €  V  we  have 

^im!A  n  z  -  it  *z|N  =  0, 

then  limj  fN(t)7iNz  -  7iNT(t)z|w  =  0  for  all  z  €  B  ,  uniformly  in  t 
N-*«>  N 

on  compact  subsets  of  [0,®). 

Proposition  2.2.  Let  {q^}  c  q  be  an  arbitrary  sequence  satisfying  qN  q*  as 

N  N  .N  N 

N  -*■  00  .  Then  the  semigroups  T  (t;q  )  and  T(t;q*)  generated  by  A  (q  ) 

and  A(q*)  respectively,  satisfy 

lim  |TN(t ;qN)z  -  T(t;q*)z|  =  0 

N-wo 

for  each  z  £  X,  uniformly  on  compact  subsets  of 


Proof.  Let  us  first  recall  some  elementary  facts  concerning  spectral 
properties  of  the  self  adjoint  operator  A(q)  (c.f.  [10,p.  291-295]  or 

[15]).  The  spectrum  of  A(q)  consists  of  a  countable  number  of  real 

00 


eigenvalues  each  of  multiplicity  1,  which  can  be  ordered 

so  that  - 00  <  ...  <  A  (q)  <  X^  ^(q)  <  ...  <  ®  and  which  are  uniformly 

OO 

bounded  above  as  q  varies  in  Q.  The  eigenfunctions  {4'^. ^ j _ 2 
A(q*)  form  a  complete  orthogonal  set  in  X.  To  apply  Proposition  2.1 
we  take  8  =  X,  8^  =  X,  ^  the  identity  operator  for  each  N  and 
A  =  A(q*),  A^  =  A^(q^).  (The  observant  reader  will  notice  that  we  could 
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equal  ly  well  take  S  =  x,  8N  =  XN,  =  P1^  with  ^  =  A(q*)  and 

=  P^A(q^),  which  in  this  case  would  be  defined  on  X^.)  Of  course, 

i)  of  Proposition  2.1  is  satisfied  trivially  while  the  stability  hypothesis 

ii)  is  the  same  as  (Hi)  which  has  already  been  established.  To  verify 
the  consistency  hypothesis  iii)  we  need  to  choose  the  set  P  .  We  define 

CO 

P  =  U spanff. .  Clearly  P  <=  dom  A(q*)  and, by  completeness  of  the 
N=1  1  n 

eigenfunctions,  V  =  8.  We  will  show  (  X-  A(q*)JP  =  P  for  any  X  in  the 
resolvent  set  of  A(q*).  Trivially  (X-A)D  c  P  holds.  To  show  that 
(X  -  A)P  =>  V  take  an  arbitrary  HP.  Without  loss  of  generality  let 
Y  be  an  eigenfunction,  say  of  A(q*).  Then  (X-A(q*))(^-y-  4\)  = 

and  therefore  Pc  (X  -  A(q*$P.  It  remains  to  show  that  A^(qN)z  -*■  A(q*)2 
for  any  z  €  P.  We  first  note  that  the  smoothness  assumed  in  (HP)  can  be 
used  to  easily  argue  that  any  (and  hence  any  z  €  P)  is  in  H  .  Thus 
for  any  z  €  P  we  have  z  €  H4  fl  hJ  . 

For  fixed  z  €  P,  let  zN  denote  P^z.  Then 

|AN(qN)z  -  A(q*)z]  =  |PNA(qN)PNz  -  A(q*)z| 

1}  £  |PN(A(qN)  -  A(q*))zN|  +  |PNA(q*)(zN  -  z) |  +  | (PN  -  I)A(q*)z| 

5  |(A(qN)  -  A(q*))zNj  +  !A(q*)(zN  -  z) |  +  | (PN  -  I)A(q*)z|  . 


For  the  second  term  we  have  the  estimate 


,  m  ,  |  q|(Dp)(x) 

l*«.*)(zN-  z)|  ixiuf0>1,  |  |K(ZN-M! 


By  Lemma  2.3  we  therefore  find  that  |A(q*)(zN  •  z)l  converges  to  zero 

like  0(N  2)  as  N  +  ,  Strong  convergence  of  PN  -*•  I  implies  that  the  third 

term  in  (2.11)  converges  to  zero.  Finally,  we  have 

,,  N.  N  qil  p,  N  p.2  \\  N  N  ,  .  , 

A(q  )z  =  — r~  (DpDz  +  pD  z  )  +  q  z  ,  which  converges  to 

*  K  Z 

Ml  2  N 

-j—  (DpDz  +  pD  z)  +  qijz,  by  Lemma  2.3  and  the  fact  that  q  -*•  q*.  This 

concludes  the  proof. 

We  have  thus  outlined  arguments  to  establish  (2.10). 


Proposition  2.3.  If  (HQ),  (HP)  and  (HF)  hold,  then  for  each  convergent 

NJ  i  N  N  i 

sequence  q‘  -*•  q*  it  follows  that  lim  |u  (t;q  )  -  u(t;q*)|  =  0  for 

N-*»  X 

each  t  €  [0,T] . 


We  further  desire  convergence  in  C(0,1)  to  obtain  the  needed 
pointwise  (in  the  spatial  variable)  convergence  to  use  with  functionals 
of  type  J.  Define  H*  =  (H*,^*,*^),  where  ,<P^=j  D<p(x)Dif  (x)p(x)dx; 

then  is  a  Hilbert  space  and  the  associated  norm  is  equivalent  to  the 

1  2  2  *^N 

H  norm  given  by  |ip  |  ,=  |Dtp(x)|  ax.  Further  we  denote  by  X  the 

0 

N  -1  N  ~N  ~N 

set  X  endowed  with  the  Hn  topology.  Let  A  (q):  X  -*•  X  be  given 

U  K 

by  A^j(q)  =  PNA(q).  Due  to  the  finite  dimensionality  of  XN  the  operator 

Aj^(q)  is  a  bounded  linear  operator  which  generates  a  semigroup  T^(t;q) 

~N  ~N 

on  X  .  To  calculate  its  exponential  bound,  we  consider  for  v  £  X 

the  following  estimate.  Letting  <*-\>2  denote  the  inner  product  in 

H^ ,  then  we  find 

<^v,  A^(q)v^>=  ^pDv,DPNA(q)v^2 

=  <pDv,  DPN(^-  D(pDv) )  +  DPNq2v]>2 
=  q2^pDv,Dvi>2  -  qi<D(pDv),PN  i  D(pDv)^2 
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=  q2^v»v^>  -  D(pL>v),  pn  ^  D(pi)v)^ 

=  q2<^v,v^  -  q^^  PND(pDv),  PN  i  D(pDv)]>x  <  q2^v,v^>  , 


since  q^  >  0.  Consequently  (see  [11 , p.85,90])  we  have 


,  N  ,  q2t,  , 

!TR(t;q)zl  1  e  I  zl  N 


for  all  z  €  X  .  This  last  estimate  implies  (in  light  of  the  weighting 


of  the  Hp  norm  by  p) 


(2.12) 


q2t  i „N 


TN(t ;q)PNz|  <  Me  | P"z|  , 
H1  H 


0 


0 


for  a  l where  M2  =  max  {p(x)|x  €  [0,1])  /min{p(x)lx  £  [0,1]). 

I  N  » 

We  next  verify  that  | P  z|  .  is  uniformly  bounded  in 

H0 

N  for  z  €  hJ  n  H2.  In  fact  |PNz|  }  =  |dPNz|2<  J D(PNz- 1Nz) j 2  «■  j»(INz- 

N  ®  3  N 

|  Uz|  2  ,  where  I'z  is  the  interpolating  cubic  spline  in  S' (A  .)  .  An 

application  of  the  Schmidt  inequality  (Lemma  2.2)  implies  that 

|PNzi  <  2  n]INz-z|7  ♦  C  N_1|d2zU  +  |Dz!,  . 

if—  *2  •-  1  *“ 


A 

with  C  a  constant.  Here  we  have  also  used  a  well  known  (see  [20,p.!>3]) 

i  N  i 

spline  interpolation  result  to  estimate  ! P ( I ‘  z-z)|7  .  Finally  for 
1  2 

z  €  dom  A(q)  =  HQ  H  II  we  obtain  (again  refer  to  [20,  p.  53]) 


(2.13) 


PN*I  !  1  X  il)2z!2  +  ]DZ  f2  ’ 


where  02  is  a  constant  independent  of  N,  We  summarize  the  above 
discussion,  in  particular  (2.12)  and  (2.13),  in 

Lemma  2.4.  If  (HQ)  and  (HP)  hold  and  z  €  H*  fl  H2,  then  for  each  T  >  0 
the  set  {T^(t;q)P^z| t  €  [0,T],  N=l,...}  is  a  bounded  subset  of  H*. 


z)l2  ♦ 


The  final  convergence  result  of  this  section  is  concerned  with  convergence 
in  the  |*|  -  norm  as  indicated  in  (2.9).  Our  results  require  severe 
restrictions  on  the  function  F  in  (2.2). 

Theorem  2.1  a)  Let  (HQ)  and  (HP)  hold  and  let  q^  q*  in  Q.  Then 

lim  ]TN(t;qN)p^--T(t;q*)zL  =  0  for  each  t  >  0  and  each  z  €  H^  0  H2 . 

w  u 

b)  Suppose  that  the  map  F  of  (2.2)  does  not  depend  on  u,  i.e.  F  =  F(q,t). 

1  2 

Suppose  further  that  t  -*•  F(q,t)  from  [0,T]  to  HQ  H  H  is  measurable 

2 

for  each  q  €  Q  and  that  q  -*■  F(q,t)  is  continuous  from  Q  to  H  for 
each  t  €  [0,T]  and  finally  for  some  k^  €  L2(0,T)  we  have  | F (q , t ) |  + 

I D2F (q, t) |  <  k, (t) ,  for  all  q  €  Q,  t  €  [0,T] .  Then  qN  >  q*  implies 
u^(t;qN)  -*■  u(t;q*)  in  C(0,1)  for  each  t  €  [0,T]  and  each  u^fq)  = 

q 4«p  e  hJ  n  h2  . 

1  2 

Proof.  To  verify  a)  let  z,y  €  HQ  n  H  .  Integrating  by  parts  and 

N  N  N  1-2 

using  the  fact  that  T(t;q*)z  and  T  (t;q  )P  z  are  in  0  H  we 
find 

!<TN(t;nN)PNz  -  T(t;q*)z,y>  J 

H0 

=  |<D(TN(t;qN)PNz  -  T(t;q*)z),  Dy)2| 

»  |<TN(t;qN)PNz  -  T(t;q*)z,D2y>2| . 

Proposition  2.2  and  Lemma  2.4  together  with  the  last  estimate  and  the 

1  N  N  N 

fact  that  dom  A(q)  is  dense  in  H^  imply  that  T  (t;q  )P  z  converges 
weakly  to  T(t;q*)z  in  H*,  uniformly  in  t  €  [0,T] .  Since  H*  is 
compactly  embedded  in  C(0,1)  by  the  Sobolev  embedding  theorem  [1  ,p.l44] 
it  follows  that  T^(t;qN)P^z-*  T(t;q*)z  in  C(0,1);  since  T  >  0  was 
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arbitrary,  this  holds  for  all  t  _>  0. 

Turning  next  to  b)  of  the  theorem,  we  note  first  that  (HF) 
is  satisfied  under  the  present  conditions  on  F.  For  each  N  we  have 

(2.14) 


N,  N.  „N,  N.  _N  ,  N, 
u  (t;q  )  =  T  ( t ; q  )P  u0(q  )  + 


N  N  N  N 

T  (t-s;qN)P  F(q  ,s)ds. 


0 


Since  |TN(t;q)PNz|  .  £  K(|d^zL  +  |Dz|  )  for  a  constant  K  independent  of 

HJ 

N=l,2,...,  q  €  Q,  and  t  £  [0,T]  by  (2.12)  and  (2. 13), the  integral  in 
(2.14)  clearly  exists  as  a  Bochner  integral  in  C(0,1).  Moreover,  for 
each  t  we  get  r(t;q>  uQ(q  )  -*•  T(t;q*)u0(q*)  and  At-SjqV  F(q  .s)  - 
T(t-s;q*)F(q*,s)  for  almost  every  s  6  [ 0  t ] ,  where  the  convergence 

is  in  the  C(0jl)  norm.  Since  the  functions  s  -*•  TN(t-s;q^)PNF(qN,s) 

1  -1/2 
from  [0,t]  to  H  are  bounded  by  the  integrable  ^unction  Kk,(*)(  min  k(x))  '  , 

U  3  xe[0,l] 

we  may  use  Lebesgue’s  bounded  convergence  theorem  when  taking  the 

limit  as  N  °°  in  (2.14)  to  find  that  uN(t;qN)  -*■  u(t;q*)  in  C(0,1) 

with  u(t;q*)  =  T(t;q*)un(q*)  +  [  T(t-s ;q*) F (q* , s)ds . 

}0 

Remark  2.1.  In  our  discussion  here  we  have  considered  only  problems 
with  Dirichlet  boundary  conditions.  Our  ideas  are  easily  extended  to 
treat  Neumann  or  mixed  boundary  conditions  as  well.  We  sketch  briefly 

some  of  the  minor  changes  required  in  the  above  preseatation.  First 

.  N 

one  certainly  must  use  different  basis  elements  (recall  how  the  B.  were 

constructed  from  the  standard  basis  elements  B^)  so  as  to  ensure 

J 

N 

X  c  dom  A(q).  This  in  turn  requires  that  one  establish  an  analogue 
of  Lemma  2.3  for  the  associated  natural  projections  P^  (through  use  of 
appropriately  defined  interpolating  splines  -  we  remind  the  reader  that 


■.t  t 
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there  are  numerous  types  of  interpolating  cubic  splines  -  see  [20, Chapter 

4] *[21]).  Finally,  minor  details  in  some  of  the  above  arguments  must 

be  modified.  For  example,  in  the  proof  of  Proposition  2.2,  the  eigenvalues 

have  multiplicity  f  2  (see  Example  4.2  of  [9]),  not  necessarily  equal  1. 

Some  of  the  integration  by  parts  arguments  (e.g.  in  establishing  the 
N 

bound  for  A^(q)  in  Proposition  2.3)  also  require  modification. 

Remark  2.2.  If  a  pure  "convection"  type  term  q^u^  also  appears  in 
the  right  side  of  (2.1),  then  one  can  use  the  theory  of  discrete  spectral 
operators  (see  Example  4.4  of  [9])  to  establish  Proposition  2.3  for  this 
case  in  much  the  same  manner  as  argued  above.  The  analogue  of  Theorem 
2.1  appears  to  be  much  more  difficult  to  obtain  however. 


-22- 


3.  Implementation  and  numerical  examples:  parabolic  equations. 

In  this  section  we  discuss  questions  related  to  the  implementation 
of  the  ideas  developed  in  section  2  and  present  some  numerical  results. 

All  of  the  results  given  below  were  obtained  using  the  fit-to-data 
criterion  J  involving  spatial  point  evaluations.  While  we  were  able 
to  establish  above  the  stronger  convergence  results  (in  C(0,1))  only  in  the 
case  of  linear  equations,  the  reader  will  see  from  the  results 
presented  below  that  the  methods  also  perform  quite  well  when  one 
uses  the  point  evaluation  criterion  J  with  nonlinear  parabolic  equations. 
The  algorithm  for  carrying  out  the  identification  of  the 

unknown  parameter  involves  two  major  tasks.  The  first  task  is: 

N  -N  N  yy 

given  N  and  u  (*;q),  find  a  q  which  minimizes  J (qju  ,y) .  This 

task  was  carried  out  by  the  standard  Levenberg-Marquardt  algorithm 

(available  in  the  IMSL  library,  routine  ZXSSQ)  and  we  shall  therefore 

not  discuss  this  part  of  the  implementation.  (For  a  discussion  of  the 

Levenberg-Marquardt  as  well  as  related  algorithms,  see  [3].) 

N  N 

The  Levenberg-Marquardt  algorithm  requires  that  J  (q)  =  J(q,u  , 

be  evaluated  for  fixed  N  at  a  sequence  of  iterates  q-1 .  Thus  we  need 

N  i 

the  approximations  u  (*;qJ)  and  computing  these  is  the  second  task. 

N 

By  (2.4)  the  values  of  u  (•;q-1)  are  obtained  by  solving  (2.6)  which 

N  N 

is  rather  easily  done  since  the  matrices  Q  and  K  appearing  there 
are  seven-banded  and  symmetric,  and  hence  they  can  be  stored  as  three 
subdiagonals  and  the  diagonal. 

Numerical  experiments  were  carried  out  for  the  general  example 
ut  =  +  q2u  +  q4f(u),  for  0  <_  x  <_  1,  and  t  >  0, 

u(0,x)  =  q3<P(x)  ,  for  0  <  x  <  1, 
u(t,0)  =  u(t,l)  =  0,  for  t  >  0, 


•<> 
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so  that  k  and  p  of  the  previous  section  are  chosen  identically  1 • 
this  case  A  =  +  q0I  and  thus  denoting  by  and  the  matrices 

with  elements 

*  <DBr  DBj>2  • 

-  «■  • 

N  N  N  N  N 

we  have  Q  =  Ay  and  K  =  -q^A^  +  q2A‘0  •  In  our  implementation,  the 

matrices  A^.A^j  were  calculated  analytically  and  stored  exactly. 

The  initial  values  and  the  nonlinear  term  require  numerical 

quadrature.  To  compute  (R^) .  =  f  ^ (x) BN(x)dx,  a  composite 

1  1  1  '0  1 

two-point  Gauss-Legendre  rule  was  employed  to  evaluate  the  integral. 

The  same  quadrature  rule  was  applied  to 

.  l  N 

(rNFL  =  ^F(q,t,uN(t;ql) ,  B^  V  =  f(q,t,x,  I  w\t;q)B^(x))B^(x)dx. 

1  1  2  JO  j=0  J  ] 

The  integration  of  the  system  of  ordinary  differential  equations 
(2.6)  was  carried  out  by  an  IMSL  routine  (DGEAR)  employing  Gear's  variable 
order,  variable  step  method.  For  the  parabolic  equations,  the  stiffly 
stable  backward  difference  methods  of  the  routine  were  used.  In  most 
parabolic  examples,  the  equations  (2.6)  are  moderately  stiff  and  while 
a  standard  Runge-Kutta  scheme  can  be  used  effectively,  it  is  more 
efficient  to  use  Gear's  stifflv  stable  method  (in  our  computations  with  this 
method,  local  error  tolerances  were  set  at  10  **) . 

N 

Finally  we  discuss  the  "inversion"  of  Q  ,  which  is  needed  not 

only  in  the  integration  of  (2.6),  (in  which  one  actually  solves  for 
N  N  N 

Q  w  (t)  even  though  w  is  used  in  the  nonlinearity),  but  also  in  the 

N  N  N 

case  where  f  =  0  to  get  u  (t;q)  from  Q  w  (t).  The  computation  of 
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N  -1 

y  =  (Q  )  z  was  carried  out  by  first  using  the  Cholesky  algorithm  to 
N  T 

decompose  Q  =  L^.L^,  ,  where  is  lower  triangular  (this  was  done  only 

once  for  a  given  N  and  the  corresponding  L^,  was  stored)  and  then 

T 

backsolving  the  equations  L^x  =  z,  L^y=x,  which  involve  only  triangular 

N 

matrices.  We  remark  that  the  banded  structure  of  Q'  is  preserved  in 

S  • 

To  demonstrate  the  feasibility  of  the  spline  approximation 
schemes  for  identification  problems,  we  generated  solutions  to  equations 
with  parameters  fixed  at  given  values  (called  "true"  values  below)  and 
then  attempted  to  identify  these  parameters  from  the  "data"  consisting 
of  values  of  the  generated  solutions . (In  some  cases  random  noise  was  added 
to  the  solution  values,  but  tests  revealed  that  this  did  not  affect 
the  efficiency  of  the  schemes.)  The  generated  numerical  solutions 
(with  the  fixed  parameter  values)  which  were  used  for  the  data  y  in 
the  fit-to-data  criterion  J  were  computed  by  an  independent  finite 
difference  method  when  closed  form  solutions  (c.g.  ir.  the  case  of  nuilinear 
examples)  were  not  readily  available. 

In  the  examples  below,  the  observations  consisted  of  the 
generated  solution  values  u(t^,Xj)  at  three  spatial  points  (x^  =  0.25, 

0.5,0.75)  sampled  at  ten  times  (t^  =  0.2,  0.4, . . . ,2.0) .  The  matrices 
C(t^,q)  in  the  fit-to-data  criterion  .1  are  taken  to  be  the  3  x  3  identity 
matrix. 


The  first  numerical  example  is  of  special  importance:  the  modal 
approximation  schemes  investigated  earlier  in  f  9  )  failed  to  perform 

satisfactorily  when  we  attempted  to  estimate  two  of  the  parameters 

(q j ,q->)  simultaneously.  The  cause  of  this  "numerical  unidentifiabi lity" 

can  be  seen  to  be  a  feature  of  the  modal  approximation  itself  and,  as 
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we  shall  see  below,  does  pot  arise  with  use  of  our  spline  approximations. 
Indeed,  in  our  numerical  investigations  we  were  unable  to  find  an  example 
in  our  class  of  parabolic  systems  for  which  the  spline  schemes  failed 
in  attempts  to  estimate  multiple  parameters. 

Example  3.1.  We  consider  the  linear  initial-boundary  value  problem 

u  =a,u  +  q„u  . 
t  41  xx  42  ’ 

u(0,x)  =  tfi(x) 
u(t ,0)  =  u(t,l)  =  0 

where  q>  (x)  =  2x  for  0  x  £  0.5  and  M>(x)  =  2(l-x)  for  0.5  £  x£  1, 

N  0 

(see  [  9  ,Ex.  6.4]).  Here  and  below  q  *  denotes  the  start-up  value 

for  the  Levenberg-Marquardt  optimization  routine,  for  the  N-th  approximate 

N 

identification  problem  (IDA  ).  Table  1  depicts  the  results  obtained 
when  q^ ,q2  are  sought  and  q^  =  1  is  assumed  to  be  known. 


TABLE  1 


N 

-N 

ql 

-N 

^2 

4 

.1020 

.8195 

6 

.1000 

.8004 

true  value 

.1 

.8 

N,0 

q 

.25 

.25 

As  we  indicated  above,  use  of  the  modal  approximations  of  [  9 ]  did  not 
produce  good  results  when  simultaneously  estimating  q^  and  q2  in 
this  example.  The  computational  findings  as  summarized  in  Table  2  reveal 
that  no  such  difficulty  is  associated  with  estimation  of  several  para- 
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meters  when  one  uses  the  cubic  spline  approximations  described  in  section 


TABLE  2 


N 

-N 

qi 

-N 

q2 

-N 

q3 

4 

.4886 

2.7856 

5.1011 

5 

.4905 

1.7779 

5.3249 

6 

.5153 

2.0594 

5.1827 

10 

.5033 

1.9884 

5.0949 

true  value 

.5 

2.0 

5.0 

N ,  0 

q 

.25 

1.0 

1.0 

Example  3.2  Next  we  consider  the  nonlinear  I-BV  problem 

Ut  =  Vxx  '  V3 

u(0,x)  =  q3ip(x), 
u(t ,0)  =  u(t,l)  =  0, 

where  cp  is  chosen  as  in  Example  3.1.  The  numerical  results  that 
were  found  are  given  in  Table  3. 


TABLE  3 


N 

-N 

^1 

-N 

q3 

4 

.4978 

5.2979 

8 

.4989 

5.1215 

16 

.4992 

5.0651 

true  value 

.5 

5.0 

qN,0 

.1 

1.0 

1.2354 

1.1055 

1.0653 

1.0 


0 
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Example  3.3.  As  our  final  parabolic  example,  we  consider  one  with  a 

different  nonlinear  equation  from  that  of  the  previous  example. 

u.  =  q,u  ♦  q  — i— 
t  M1  xx  M4  1+u 

u(0,x)  =  q3<P(x), 

u(t ,0)  =  u(t,l)  =  0  , 

where  q>  is  chosen  as  in  the  previous  examples.  The  numerical  findings 
are  recorded  in  Table  4. 


TABLE  4 


N 

-N 

<*1 

-N 

^3 

-N 

*4 

4 

.  5309 

5.2963 

1.0781 

8 

.5126 

5.1208 

1.0316 

16 

.5067 

5.0648 

1.0170 

true  value 

.5 

5.0 

1.0 

N,0 

q 

.1 

1.0 

0 
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4.  A  class  of  hyperbolic  partial  differential  equations. 

In  this  section  we  consider  identification  problems  for  the 
one  dimensional  nonlinear  hyperbolic  equation 

(4.1)  u  =  q1D2u  +  q2ut  +  q3u  +  f (q6,t,x,u,ut) ,  for  t  >  0,  x  €  [0,1], 

with  boundary  and  initial  conditions, 

(BC)  u(t , 0)  =  u (t , 1 )  =  0,  for  t  >  0 
(IC)  u(0,x)  =  q4<P(x), 

ut(0,x)  =  qs^(x),  for  x  €  [0,1]. 


We  could  allow  equally  well  for  multiple  unknown  parameters  in  the 


initial  conditions,  as  for  example  u(0,x)  =  .1  q*4>.  (x)  . 

Following  the  procedure  of  section  2,  we  rewrite  (4 . 1) - (BC) - (IC) 


as  an  abstract  evolution  equation  in  an  appropriate  Hilbert  space. 

'2( 


2  0  2 
Standard  results  imply  that  D  in  H  =  L^(0,1;R)  with  dom  (D  ) 


1  2  S  2  2 

Hg  D  H  is  a  selfadioint  operator  satisfying  ^-D  z,z^  —  lz  I? 

2 

for  every  z  £  dom(D  ).  We  again  assume  that  (HQ)  holds,  except  now 

modified  in  that  Q  c  R6  .  With  (HQ)  holding  we  define  V(q^)  =  (Ho’^*’*^V(q  p 

where  we  endow  H?.  with  the  topology  defined  bv  the  inner  product 
u  -1 


<w,z> 


V(qt) 


q^Dw(x)Dz(x)dx  =  ^q^Dw.Dz^- 


Clearly  V(q^)  is  a  Hilbert  space.  Define  X(q)  =  V (q j )  *  H°  with  the 

usual  product  topology  generated  by  *  ^z\  ,z2^=  ^W1  ’  21^  V (q^) +  ^2*  Z2^2 

Once  again,  </,•>  and  |.|  (or<<,/>q  and  |.[q  whenever  dependence  on 
the  vector  q  must  be  emphasized)  denote  the  inner  product  and  associated  norm  ! 
X  throughout  this  section.  We  may  now  rewrite  (4.1)  in  X  (taking 
v  =  ut)  as 
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(•[*])  "  A(q)(y[^)+  F(q,t,u(t)  ,v(t)) ,  for  t  >  0, 

(4.2) 

/u(0)\  =/q44>\ 

\v(0)/>  V  4>/  ’ 

where  (tp ,^)  £  X,  dom  A(q)  =  hJ  fl  H‘  x  H'  ,  A(q)  =^q  D2+<^ 

and  F  is  the  operator  given  by  F(q,t ,u(t) ,v(t) )  =J  0  \ 

yf  (q»t, • ,u(t , • ) »v(t , »))y 

As  in  section  3  we  consider  the  identification  problem 

A  (y 

(IDA)  Minimize  J (q,u(* ;q) ,y)  over  q  £  Q  c  R  subject  to  u(*;q) 
satisfying  for  t  £  [0,T] 


(4.3) 


/u(t;q)\  /u(0)\  +  [t 

^v(t;qy  Ht,qJlv(0)l  J0 


T(t-s;q)F(q,s,u(s;q)  ,v(s;q))ds. 


Here  T(t;q)  denotes  the  semigroup  in  X(q)  generated  by  A(q) 

(see  [9]).  We  shall  assume  that  (HF)  with  u  replaced  by  (u,v) 
holds  throughout  this  section;  this  suffices  to  ensure  existence  of 
a  unique  solution  of  (4.3)  on  [0,T] .  Again,  conditions  on  f  (q^.  ,t  ,x,u,ut) 
that  imply  (HF)  can  be  found  in  [9].  Since  the  first  component  u(t) 
of  the  solution  vector  (u(t),v(t))  lies  in  H*,  point  evaluations  of 
u(t)  clearly  pose  no  difficulties  in  this  example.  Moreover,  a  para¬ 
meter  dependent  convergence  result  in  the  norm  of  the  state  space 

analogous  to  (2.10)  will,  in  this  case,  be  sufficient  to  prove  that 
-N 

solutions  q  of  the  approximating  identification  problems  and  the  associated 
N  -N 

fit-to-data  term  J  (q  )  converge,  respectively,  to  a  solution  q  of  the 
original  problem  (IDA)  and  its  associated  term  J(q). 

An  interesting  aspect  of  this  example  is  that  now  the  norm 
of  the  state  space  depends  on  the  parameter  vector  q.  The  desirability 
for  allowing  for  such  generality  results  from  the  fact  that  the 
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operator  A(q)  will  be  dissipative  in  this  weighted  inner  product.  In 
fact,  with  (HQ)  obtaining,  one  can  easily  establish  the  existence  of 
a  constant  w  independent  of  q  G  Q  such  that 

(4.4)  <^A(q)z,^  <  u|  z|q 

for  all  z  €  X  (See  [  9 ]  .  One  can  also  rescale  the  state  variable 

to  avoid  the  state  space  dependence  on  q,  but  this  does  not  lead  to 

any  simplification  in  the  computational  scheme  obtained.) 

We  now  define  the  state  approximation  scheme.  As  in  the 
3  N  N  N 

previous  example,  define  Sp(A  )  =  span  {Bq,...,B^  }.  For  the 

approximating  subspaces  we  take  XN  =  S^(AN)  *  S^(AN)  cz  H*  n  *  H*  = 

N  N  N 

dom  (A  (q) ) .  In  terms  of  basis  elements,  X  (q)  =  span  , . . . >32n+2^  ’ 
where 

(B^1_1,0)T  ,  1  <  i  <  N+l 
(°'Bi_N_2)T  '  N+2i  1  £  2N+2  * 


We  again  seek  an  approximation  to  the  state  of  the  form 

(uN(t;q),  v^ft ;q) )T  =  £  w^(t)8^  , 

i=l 

by  using  a  Galerkin  approach.  Proceeding  as  before  we  arrive  at  a 
semidiscrete  approximation  to  (4.2)  by  a  system  of  2N+2  ordinary 
differential  equations: 


QNwN(t)  =  KNwN(t)  ♦  RNF(q,t,uN(t),vN(t)) 

v  V  M  T 

QV(0)  =  RN((q4d>,q^)‘), 

where  wN(t)  =  col  (w^ (t)  , . . .  ."^+2^))  »  and 
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A.,  <v  *?>, 

Aj  =  <e".  A«t)ej>q 

(A);  =  <e“.  *>q  . 

we  again  find  it  more  convenient  to  write  this  in  terras  of  the 

orthogonal  projections  P*^ (q) : X (qj)  -*■  X^(q),  and  defining  A^(q)  =  P^CqjMqJP  (q)  > 
we  obtain 

wN(t)  =  AN(q)wN’(t)  +  [PNF(q,t,uN(t)  ,vN(t)  j] 

(4-5)  M  M  T 

wN(0)  =  [PN(q4<P,q5<JO  *], 

where  AN(q)  is  the  matrix  representation  of  AN(q)  and  [• ]  denotes 

N 

the  coordinate  vector  of  an  element  in  X  .  Note  here  that  as  a  con¬ 
sequence  of  the  parameter  dependent  norm  the  projections  also  depend 
on  q. 

As  before  we  solve  (IDA)  iteratively  by  solving  the  sequence 

of  approximating  problems  (IDA^)  with  (2.8)  replaced  by  (4.5).  Using 

considerations  similar  to  those  of  section  2  it  is  easily  seen  that 

solutions  qN  of  the  approximating  problems  exist  and  that  any  limit 

q  of  a  convergent  subsequence  of  q^  is  a  solution  of  (IDA) ,  provided 

N 

that  for  any  arbitrary  sequence  q  in  Q  we  have 
(4.6)  lim  |qN-q*|  =  0  implies  lim  | (uN (t ;qN) ,vN (t ;qN) )  -  (u(t ;q*) , v(t ;q*) ) |  =  0, 

ty-joo  A 

for  each  t  €  f 0 ,T] .  Employing  the  mild  forms  of  equations  (4.2)  and 
(4.5)  satisfied  by  (u(t ;q*) ,v(t ;q*))  and  (uN(t ;qN) ,vN(t;qN))  respectively, 
one  can  show  as  before  that  (4.6)  holds  provided  that  conditions  (Hi)-(Hiii) 
explained  in  section  2  can  be  verified.  Since  (Hi)  will  be  discussed 
together  with  (Hiii),  we  immediately  turn  to  (Hii).  In  the  calculations 
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below  we  will  freely  use  the  fact  that  under  (HQ)  all  the  X(q)  norms 

are  equivalent  as  q  varies  in  Q. 

N 

To  show  that  P  -*•  I  strongly  in  X,  we  will  need  to 
modify  the  previous  spline  approximation  Lemma  2.3  to  obtain  V(q^)  norm 
estimates.  First  we  state  a  spline  interpolation  estimate  in  the  V(q^) 
norm;  (see  Lemma  2.1,  of  which  this  is  simply  a  restatement  of  the  last 
two  inequalities) . 


Lemma  4 . 1  Let  INz  be  the  cubic  spline  interpolant  (from  S3(AN))  to 


z  €  H  .  Then 


tN 

I  z  -  z 


I  D(I  z  -  z) 


,,)  -  <iJ/2c1n's|d4z|2  , 
•(q,)  ^  qJ/2C2N-2(D4z  |  2  - 


Next  we  establish  a  V(q^)  norm  Schmidt  inequality  for  splines. 


Lemma  4.2  Let  s€S3(AN).  Then  |s|v^  ^  N  q|/z|  s|  2 

U’Cqj)  -  ^2  Nq/  2  =  NlSlv(q1)  ' 


1/2. 


and  Ds 


Proof.  Take  s€S3(AN).  Then 


1  N  ,x.  N  rx- 

|Ds(x)rdx  =  l  |Ds(x)i  dx<  I  k.N  |s(x)|  dx 

0  i=l  'x.  .  i  =  l  A  '  x.  . 

i-l  i-l 


where  we  have  used  Lemma  2.2.  This  yields  the  first  inequality.  Since 
Ds  is  a  polynomial  of  degree  2  on  [x^  ^ , jc^.  ] ,  a  similar  argument  provides  the 
second  estimate. 

Finally  we  have  the  desired  projection  estimates  in  the  V(qp 


norm. 
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N  1  3  N 

Lemma  4.3.  Let  S^(A  )  be  the  orthogonal  projection  in  the 

4  1 

V(qj)norm.  Then  for  z  €  H  fi  we  have 

|PjZ  -  z|v  )  <  qJ/2C1N'3|D4z|2  ,  and 


D(Pr  '  Z)U(qj)  <  qJ/2(2^C1  +  C2)N'2|D4z|2  . 


Proof.  The  first  estimate  follows  directly  from 


I”?1  -  *!»«,,)  i  I1"2  -  z'vtqi3  i  qJ/2c1N-3|D4zi2  , 


N  3  N  3  N 

since  the  interpolate  I  z  in  S^(A  )  is  the  same  as  that  in  S  (A  ) 

4  1 

whenever  z  €  H  (1  -  see  Lemma  4.1.  To  obtain  the  second  estimate, 

we  first  write 

lD(^2  -  2>lv(qi)  ilD<=  -  iN^Ivcci3  *  I011"2  -  . 

The  desired  estimate  on  the  second  term  follows  directly  from 

the  second  estimate  in  Lemma  4.1.  To  estimate  the  first  term,  we  use 

the  Schmidt  inequality  (from  Lemma  4.2)  along  with  Lemma  4.1  to  obtain 

<  ^N[|P^  -  z|v(qp  *  P  -  I^lv(qi,] 

5  2^N|INz  -  i|v(  .  5  2PT'N'2C,q|/2|D4l|2  . 

N 

Proposition  4.1.  Assuming  (HQ),  the  projections  P  converge  strongly 
to  I  in  X. 

Proof.  Let  z  =  (zrz2)  €  H4  n  hJ  *  H4  fl  hJ.  Then  PNz  =  (z^,z2)  s  (P^  .P^) , 
where  Pj  is  the  projection  of  the  first  coordinate  into  Sq(A^)  in 


4 
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N 

the  V C M j )  norm  and  is  the  projection  of  the  second  coordinate  into 

Sq(AN)  in  the  H®  norm.  Then 

l dN  ,  2  ,  nN  .2  , nN  ,2 

lP  2  ~  z!  ■  lPl2l  '  2ll  VCq,)  *  |P2Z2  -  z2>2 

<  (q{/2c1N-3|D4t1|2)2  ♦  (C0n"4|d4z2|2j2  . 

N 

This  estimate  together  with  the  boundedness  of  P  imply  that 
PNz  ■+  z  for  any  z  €  X. 

Finally  we  discuss  the  parameter  dependent  convergence  of 
the  linear  semigroups. 


Theorem  4.1. 


Let  (HQ)  hold. 


Then  the  semigroups  T(t;q)  and 


N  N 

T  (t;q)  generated  by  A(q)  and  A  (q)  satisfy  exponential  bounds 


|| T (t  ;q)  ||  _<  eWt  and  1 1 T1^ (t ; q)  ||  £  for  some  real  u  independent 

of  N  and  q.  Moreover  for  any  sequence  {qN}  converging  to  q*  in 
Q  we  have  |TN(t;qNjz  -  T(t;q*)z|  ■>  0  uniformly  on  [0,T]  for  each 


wt 


e  x. 


ii  II  U)t 

Proof.  The  bound  ||T(t;q)|]  £  e  is,  of  course,  a  consequence  of 

N  N 

(4.4).  As  in  section  2,  A  (q)  clearly  generates  a  semigroup  T'  (t;q) 

whose  exponential  bound  is  seen  directly  from 

<AN(q)z,z>q  <  <'A(q)PN(q)z,  PN(q)z>q  <  w|PN  (q  )  ;  I  2  <  w  j  2 '  ‘  . 

N  N 

Here  we  used  the  fact  that  P  (q):X(q)  X  (q)  is  the  orthogonal  proiec 
t  ion . 

We  now  turn  to  Proposition  2.1  to  establish  the  convergence 
result  of  this  theorem.  Let  S  *  X(q*),  B^=  XN(qN)  and  A  =  A(q*),  A^  »  A^*(qN 
Then  ii)  of  Proposition  2.1  holds  for  the  family  of  semigroups  T^(t;q^). 
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Letting  ti^ : X (q* )  -*■  X(q^)  be  the  canonical  isomorphism  between  X(q*) 

and  X^q^),  it  follows  immediately  that  |i^z|  ■+  |z  |  from  the  hypothesis 

that  q1^  ■>  q*  and  thus  i)  of  Proposition  2.1  is  also  satisfied. 

To  define  V  >  we  note  that  <J>j  (x)  =  (/2? ) sin  (px)  and 

$.  (x)  = /T  sin  (jnx),  j  =  1,2,...,  form  complete  orthonarmal  sets  for 

V(l)  and  respectively.  Let  V^(x)  =  ($j,0)T,  for  j  =  1,...,N,  and 

Yj  =  (0.*  )*  for  j  =  N+l , . . .  ,2N.  Then  forms  a  complete 

orthonormal  set  for  X(q),  where  q  =  (1,0,...,0),  and  a  complete 

orthogonal  set  for  X(q)  for  q  arbitrary.  We  take  =  span{Y^, . . .  Y ^  } 

and  V  =  U^.  For  this  choice  of  p  it  clearlv  follows  that  V  c  dom  A(q*) 
N  =  1 

and  V  =  B  •  Also,  from  the  definition  of  A(q)  it  is  easily  argued 
that  for  X  >0  sufficiently  lar^e  (Xl  -A(q*))P  =  V  ,  so  that 
(AI  -  A(q*))P  is  dense  in  8. 

Finally,  to  establish  the  convergence  part  of  the  consistency 

N 

hypothesis  iii)  in  Proposition  2.1  (suppressing  the  notation  tt 

for  the  canonical  isomorphism),  we  see  that  for  each  z  =  (z^,z2)€  V  . 

once  again  (see(2.11)) 


(4.7) 


|AN(qN)z  -  A(q*)z|  =  ] PN(qN)A(qN)PN(qN)z-  A(q*)z| 

<  !(A(qN)  -  A(q*))PN(qN)z|  +  | A(q*) (PN (qN) z  -  2) |  +  | (PN (qN) - 1 )A (q*) 1 1 • 


These  three  terms  will  now  be  estimated  separately.  First  we  write 
explicitly  the  second  term  as 

(N  \ 

Z2  2  \ 

q*D2(Zj  -  zx)  >  q^Cz^-Zj)  +  <\*{ZN2-Z2)  I  . 


,N,  N.  . 
(P  (q  )z-z) 


-36- 


Thus 

jA(q*)  (PN(qN)z  -  z)\2  =  i  z‘2  -  z2iv(q|)  +  '  2P  +  ^3^  “  2P  +  ^2(z2  “  Z2^ 

<.q*|D(Z2  -  Z2)  i  2  +  ^i!d2(Zj  -  Z1)f2  +  |  q*  i  i  2^  -  Zj  (  2  +  IqjUz^  “  Z2^2^  ' 

Using  the  fact  that  *^cjj  | ( z‘^  -  zpi2  =  l^zi  '  zl^V(q*)  ’  we  see  t^ia<: 
this  last  estimate  together  with  Lemma  2.3  and  Lemma  4.3  imply  that 
|A(q*)  (P^(q^)z  -  z)]  -*•  0  as  N  -*•  00  . 

N  N 

For  the  last  term  in  (4.7)  we  have  convergence  from  P  (q  )  -*•  I 
strongly  on  X.  Finally  convergence  of  the  first  term  on  the  right  in 
(4.7)  can  be  seen  from 

/  0 

(A(qN)  -  A(q*))PN(qN)z  =  {  .  N  „*..n2  N  ,  N  ..  N  ,  N  ..  N 
M  M  M  y(qx-qpD  Z1  +  (q3-q*)z1  +  (q2-q*)z2 


,2N 


N 


N 


and  the  fact  that  I)  z^  -*•  D  z^,  z^,  and  z7  -*■  z2  in  H  and  q  •*  q’ 

Thus  Proposition  2.1  is  applicable  to  establish  the  desired 
convergence  in  the  theorem. 
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5.  Implementation  and  Examples:  the  Hyperbolic  Case. 

This  section  is  devoted  to  a  discussion  of  the  computer 

implementation  of  the  theory  developed  in  the  previous  section  and  to 

the  documentation  of  some  of  our  numerical  findings.  As  in  the  parabolic 

case,  the  optimization  part  of  the  algorithm  was  carried  out  using  the 

IMSL  implementation  of  the  Levenberg-Marquardt  algorithm.  For  given 

values  of  N  and  q  the  state  approximations  must  be  obtained  by 

solving  (4.5).  Again  the  matrices  involved  in  this  equation  have  a 

structure  which  permits  efficient  solution  of  the  system  of  the 

approximating  ordinary  differential  equations.  The  required  matrices 
N  w  _  1  N 

for  A  (q)  =  (Q‘  )  K  are  given  by 


where  (Q1^) .  .  =  (K1?).  .  =  <"B‘N,BN  .  ,  (Q^) .  .  =  /bN,BN/>  > 

^ri,j  v  1  i,l  N  i  j/^V(q1)’  2  l ,1  >  l  j'2 

(K*).^  =  ^lD2B?+q5Bi’Bj>2  ’  and  K3  =  q3Q2  *  Equ^alently, 
if  we  let  (.^)j  DB^^^and  (A^).  ,  these  matrices 


become 


and 


aN 

<tiai 


0 

0 


kn  =1 


i  aN  an  an 
\qlVq3A2  q2A2 


The  matrices 


aN  . 
A  and 


are  exactly  those  which  were  used  in  the 
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parabolic  case  and  hence  they  were  computed  analytically  and  stored. 

N  N  - 1  N  . 

The  "inversion"  of  Q  in  computing  (Q  )  K  w  is  carried  out  as 


follows:  We  have 


(QN)_1KN  = 


.N.-l 

«1  (A1} 


(A*)'1 


0 


an 

qlAl 


an  .n  .n 

-qiAl+q3A2  q2A2 


^iCAj'^l  q2I 

Given  w  =  col(w  ,w2)  where  »1>V2  are  vectors  of  length  N+l,  we 

~  ~  N  -1  N  N 

compute  w  =  col(w  ,w2)  =  (Q  )  K  w  by  first  computing  y  =  A^  . 

Then,  since  is  stored  as  its  Cholesky  factor  L^,  where 

a2  =  LnlJ,  we  compute  by  solving  LNz  =  y,  L^r  =  z, 

for  r.  Finally  we  have  w  =  col (w2,-qjr+q3Wj+q2w2) . 

The  projection  of  the  initial  data  and  the  nonlinear  term 

again  require  numerical  quadrature.  Recall  that  the  first  coordinate 

of  the  projection  map  is  projection  with  respect  to  the  V(q^)  inner 

product.  Given  the  initial  data  u(0,x)  =  q^p(x),  v(0,x)  =  q^Cx), 

we  compute 

qA  ( (Q1?)”1  0  \  (vC*’B?>V(q  ) 

CQ  "A'  1 


(QN)-1RNI 


(Q^)"1 
o 


d2<>2 

where  <^.d2b^>2  and  are  computed  via  a  composite 


N  N 

Gauss -Legendre  rule  and  the  "inversions"  of  ,  An  are  carried  out 

through  the  Cholesky  decomposition.  The  projection  of  the  nonlinear 

term  is  computed  in  the  same  manner.  We  point  out  here  that  for  this 

N  N 

class  of  hyperbolic  equations  the  projection  operator  P  (q):X(q)  -*■  X  (q) 
actually  does  not  depend  on  q,  as  can  be  seen  from  the  above  cal¬ 
culations  . 

The  approximating  ordinary  differential  equations  were  again 
solved  using  the  IMSL  package  DGEAR,  but  for  the  hyperbolic  examples 
the  Adams  multistep  option  of  that  package  was  chosen. 

"Data"  were  generated  as  in  the  parabolic  case  above  (using  an 
independent  finite  difference  scheme  to  generate  solutions  corresponding 
to  "true"  parameter  values).  Except  where  otherwise  noted,  the  obser¬ 
vations  used  in  the  following  examples  consisted  of  the  solution  values 

■A 

u(t.,Xj)  at  three  spatial  points  (x^  =  .25,. 5,. 75)  sampled  at  ten  times 
(t.  =  .2, .4,... ,2.0). 


Example  5.1.  The  first  example  for  hyperbolic  equations  involves  the 
linear  ini tial -boundary  value  problem 

utt  =  Vxx  +  Vt  ’ 

u(t,0)  =  u(t,l)  =  0  , 
u(0,x)  =  2x(l-x), 
ut(0,x)  =  0  . 

The  numerical  results  can  be  found  in  Table  5,  where  the  same  notation 


as  in  section  3  is  used. 
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TABLE  5 


N 

-N 

<*2 

4 

1 .9967 

-.9984 

true  value 

2.0 

-1.0 

N ,  0 

q 

1.4 

.4 

Our  numerical  results  are  given  in  Table  7.  Note  that  the  estimates 
are  very  good  at  N=8,  while  numerical  error  in  implementing  the  scheme 
(e.g.  inverting  large  matrices)  is  apparent  at  N=16. 
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TABLE  7 


-N 

-N 

N 

ql 

% 

4 

1.9894 

2.0538 

8 

1.9992 

.9410 

16 

1.9979 

1.0833 

true  value 

2.0 

1.0 

N ,  0 

q 

1.0 

.5 

Example  5.4.  Here  we  consider  the  nonlinear  problem 

Utt  ~  Vxx  Ut  +  q6  1+u  ’ 
u(t,0)  =  u(t,l)  =  0 

u(0,x)  =  2x(l-x), 

ut(0,x)  =  0. 


The  observations  for  the  results  detailed  in  Table  8  consisted  of 
5  spatial  points  and  10  time  samples,  as  in  the  previous  examples. 

This  example  was  also  run  using  only  two  time  samples  (t  =  .2,. 4) 

at  the  3  spatial  points  and  rapid  convergence  of  the  estimated  parameters 

to  the  true  values  was  also  observed  in  this  case. 


TABLE  8 


N 


-N 

*1 


-N 


4  2.0084 

8  2.0003 

true  value  2.0 

CN’°  1.0 


3 . 0866 
3.0128 
3.0 
1.0 


Example  5,5.  A  third  type  of  nonlinearity  was  used  in  the  example 

utt  =  Vxx  +  Vt  +  %  sin  11 

u(t,0)  =  u(t,l)  =  0, 
u(0,x)  =  2x(l-x), 


ut(0,x)  =  0. 
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Tht'  reader  will  find  the  results  in  Table  9. 


TABLE  9 


N 

-N 

ql 

-N 

^2 

-N 

^6 

4 

1.9549 

-.9994 

2.5559 

8 

2.0023 

-.9992 

3.0351 

16 

1.9984 

-.9998 

2.9969 

true  value 

2.0 

-1.0 

3.0 

N ,  0 

q 

1.0 

0 

1.0 

Example  5.6.  As  a  final  example  we  present  in  Table  10  the  results  for 

Utt  =  Vxx  +  ^3U  +  %  ’ 

u(t,0)  =  u(t,l)  =  0 
u(0,x)  =  2x(l-x) 
ut(0,x)  =  0. 

TABLE  10 


N 

-N 

**1 

-N 

^3 

-N 

% 

\ 

1.9653 

.6509 

2.9892 

8 

2 . 0055 

1.0516 

2.9988 

16 

1.9999 

.9971 

3.0013 

true  value 

2.0 

1.0 

3.0 

N  ,0 

q 

1.0 

0 

1.0 
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6.  Concluding  Remarks. 

In  this  paper  we  have  shown  that  spline  approximations 
may  be  profitably  used  to  develop  schemes  for  estimation  of  unknown 
parameters  in  initial-boundary  value  problems  for  second  order  partial 
differential  equations.  The  practical  utility  of  our  ideas  is  supported 
by  our  computational  experience  with  a  large  number  of  examples,  a 
summary  of  which  is  also  given  in  this  paper.  The  use  of  spline  functions 
in  the  context  of  parameter  identification  problems  is  not  new;  see,  for 
instance,  references  found  in  the  survey  articles  [2],  [12] , [16] , [18] . 

It  is,  however,  our  belief  that  ours  is  the  first  presentation  of  a 
complete  theoretical  treatment  (i.e.  convergence  proofs  for  the  parameters, 
optimal  states  and  optimal  fit-to-data  values)  for  spline-based  methods 
for  a  large  class  of  equations  along  with  reports  on  a  careful  numerical 
testing  of  the  methods  on  examples. 

The  fundamental  ideas  -  involving  use  of  a  semigroup  theoretic 
framework  for  the  approximation  of  identification  and  control  problems 
governed  by  partial  differential  equations  -  which  are  the  basis  of  the 
convergence  results  of  the  present  paper  were  first  announced  in  [8]. 

A  complete  presentation  of  details  of  the  approximation  framework  was 
later  given  in  [9].  In  fact  the  present  paper  is  a  companion  paper  to 
[9]  detailing  our  work  on  spline  approximations  that  was  performed 
simultaneously  with  our  efforts  on  modal  approximations;  a  report  on 
the  latter  is  contained  in  [9]  as  an  application  of  the  general  framework 
developed  there. 
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The  essential  step  in  our  considerations  above  and  in  [9]  in 

approximation  of  the  abstract  differential  equation  is  the  approximation 

N  N 

of  the  infinitesimal  generator  A  by  P  AP  .  The  importance  of  this 
classical  approximation  of  an  unbounded  operator  by  a  sequence  of  bounded 
operators  to  problems  in  system  theory  using  semigroup  methods  (in  par¬ 
ticular  the  Trotter-Kato  theorem  -  Prop.  2.1  above)  was  first,  we  believe, 

N  N 

pointed  out  in  [7],  where  the  P  AP  scheme  was  used  with  spline 
approximations  for  functional  differential  equations. 

In  this  paper  we  chose  two  specific  classes  of  examples  to 
illustrate  use  of  our  spline  approximation  ideas.  Particular  boundary 
conditions  and  fit-to-data  criteria  were  selected.  However,  given  the 
general  results  of  [9]  and  the  technical  developments  discussed  above, 
it  should  be  obvious  that  many  rather  easy  but  important  generalizations 
are  possible.  We  mention  a  few  of  these  which  we  have  already  investigated 
in  applying  our  ideas  to  a  number  of  specific  problems.  First,  methods 
for  problems  with  mixed  boundary  conditions  involving  unknown  parameters 
can  be  readily  obtained  -  see  [5]  for  a  partial  discussion  of  some  of 
our  results  in  this  direction.  Furthermore,  we  have  previously  discussed 
elsewhere  [4]  how  one  might  use  our  ideas  to  develop  techniques  for 
estimation  of  unknown  variable  coefficients  in  partial  differential 
equations.  Indeed  we  have  already  successfully  used  spline  methods 
along  the  lines  developed  above  to  estimate  spatially  varying 
transport  coefficients  in  insect  dispersal  models  from  field  data  (in 
collaboration  with  P.  Kareiva) .  (These  efforts  are  detailed  in 
separate  manuscripts  currently  in  preparation.)  In  related  efforts 


a  rather  different  approach  (not  involving  linear  semigroup  approxination 
results)  is  taken  in  [6]  where  spline  methods  for  estimation  of  time 
varying  coefficients  are  developed. 

The  framework  of  [9]  and  the  present  paper  can  also  be  used 
to  treat  higher  order  equations  such  as  those  arising  in  the  study  of 
large  space  structures  and  other  elastic  bodies  (see  [4],  [5]).  In 
particular  we  have  in  this  spirit  developed  methods  employing  cubic 
and  quintic  spline  approximations  for  equations  such  as  those  arising 
in  both  the  Tuler-Bernoulli  and  Timoshenko  theories  for  beams.  The 
ideas  and  some  of  our  numerical  findings  are  reported  in  [5] . 

In  a  modification  of  the  theory  presented  here  (not  a  PNAPN 
type  scheme)  but  using  essentially  the  same  semigroup  approximation 
framework,  we  have  developed  methods  that  permit  one  to  use  lower  order 
splines  (e.g.,  piecewise  linear  splines  for  second  order  equations)  in 
place  of  the  cubic  splines  employed  in  this  paper.  A  preliminary 
discussion  of  theoretical  and  computational  results  for  second  order 
equations  can  be  found  in  [13].  A  complete  discussion  of  these  ideas 
for  second  and  higher  order  equations  along  with  numerical  results 
will  be  published  in  a  manuscript  currently  in  preparation. 

All  of  the  specific  results  in  this  paper  and  those  mentioned 
above  to  appear  elsewhere  entail  spline  approximations  for  problems 
involving  one  dimensional  spatial  domains  (x  £  R*).  However,  we  have 
investigated  extensions  of  our  ideas  to  treat  problems  in  higher 
dimensional  domains  (x  €  Rn)  and,  not  surprisingly,  there  appears  to  be 
little  difficulty  in  making  the  needed  extensions  of  the  theory  (e.g., 
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see  some  of  the  estimates  for  approximation  using  product  basis  elements 
in  Chapter  6  of  [20]).  We  are  currently  considering  the  computational 
feasibility  of  the  spline  methods  for  such  problems,  but  expect  once 
again  to  find  that  they  perform  well. 

Finally  we  note  that  the  results  in  this  paper  and  the 
framework  of  [9]  concern  development  of  schemes  that  involve  semi¬ 
discretization  in  that  one  approximates  the  original  partial  differential 
equation  by  a  sequence  of  ordinary  differential  equations.  One  can  also 
develop  a  theoretical  framework  to  combine  spline  approximations  in 
the  spatial  coordinates  with  time  discretizations  resulting  in 
approximation  by  a  sequence  of  difference  equations.  In  particular 
the  theory  developed  in  [19]  involving  full  discretization  methods  for 
functional  differential  equations  can  be  used  to  develop  analogous 
techniques  for  partial  differential  equations.  We  are  currently  pursuing 
investigations  of  these  ideas. 
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